mzM 


AD-A093  299  MISSISSIPPI  STATE  UNIV  MISSISSIPPI  STATE  ENGINEERING— ETC  F/G  12/1 

NUMERICAL  GENERATION  OF  TVO-DIMENSIONAL  ORTHOGONAL  CURVILINEAR  — ETCtU) 
JUN  80  2  U  VARS I  *  R  A  VEED»  J  F  THOMPSON  AFOSR-8O-O105 

UNCLASSIFIED  MSSU-EIRS-ASE-BO-3  AFOSR-TB-AO-1 33Q  mi 


AFOSR-YR-  80*1339 

NUMERICAL  GENERATION 


TWO-DIMENSIONAL  ORTHOGONAL 
CURVILINEAR  COORDINATES 

I 

IN  AN 

EUCLIDEAN  SPACE  * 


% 


ENGINEERING 


INDUSTRIAL  RESEARCH  STATION 


Department  of  Aerospace  Engineering 


Z.  U.  A.  Warsl 


R.  A.  Weed 


Thompsc 


DTIC 

ELECTE 
OEC  30  1980 


Mississippi  Stats  University 
Mississippi  8tats,  Mississippi  39762 


jvn  os  mo 

MSSU-EIRS-ASE-80-3 


COLLEGE  OF  ENGINEERING  ADMINISTRATION 


WILLIE  L.  MCDANIEL,  PH.D. 

DEAN.  COLLEGE  OF  ENGINEERING 

WALTER  R.  CARNES,  PH.D. 

ASSOCIATE  DEAN 

RALPH  E.  POWE,  PH  D. 

ASSOCIATE  DEAN 

LAWRENCE  J.  HILL,  M.S. 

DIRECTOR  ENGINEERING  EXTENSION 

CHARLES  0.  CLIETT,  M.S. 

AEROSPACE  ENGINEERING 

WILLIAM  R.  FOX,  PH.D. 

AGRICULTURAL  &  BIOLOGICAL  ENGINEERING 

JOHN  L  WEEKS,  JR.,  PH  D. 

CHEMICAL  ENGINEERING 

ROBERT  M.  SCHOLTES.  PH.D. 

CIVIL  ENGINEERING 

B.  J.  BALL,  PH.D. 

ELECTRICAL  ENGINEERING 

W.  H.  EUBANKS,  M  ED. 

ENGINEERING  GRAPHICS 

FRANK  E.  COTTON,  JR.,  PH.D. 

INDUSTRIAL  ENGINEERING 

C.  T.  CARLEV,  PH.D. 

MECHANICAL  ENGINEERING 

JOHN  I.  PAULK,  PH.D. 

NUCLEAR  ENGINEERING 

ELDRED  W.  HOUGH,  PH  D. 

PETROLEUM  ENGINEERING 


For  additional  copies  or  information 
address  correspondence  to 

ENGINEERING  AND  INDUSTRIAL  RESEARCH  STATION 
DRAWER  DE 

MISSISSIPPI  STATE  UNIVERSITY 
MISSISSIPPI  STATE  MISSISSIPPI  39762 

TELEPHONE  1 601 )  325  2266 


Mississippi  State  University  does  not  discriminate  on  the  basis  of  race  coioi  religion  national  onqm 
sex.  age.  or  handicap 

In  conformity  with  Title  IX  of  the  Education  Amendments  of  >972  and  Section  504  of  the  Rehabilitation 
Act  of  1973.  Dr  T  K  Martin  Vice  President,  610  Allen  Hail  P  O  Drawer  j.  Mississippi  Stale  Mississippi 
39762  office  telephone  number  325-3221  has  been  designated  as  the  responsible  employee  to 
coordinate  efforts  to  carry  out  responsibilities  and  make  investigation  of  complaints  relating  to 
nondiscrimination 


SECURi 


PORT  DOCUMENTATION  PAGE 

n 


llgnnnT  TTMi.i'rirr  _ ,  y_ _  >  12.  govt  accession  no.i  3  •  mcui^icn  »  ar-  «  1  p  i\»  wwiiippii 

AF0SRlSTRl8 (f -  1  3  39; «  0  a  r1  fj  # 


b. 


i 

J 


I  WA  uAULinfluif.  Af.£urv.Ji4uF  A  annRf  SVif  jfiUayn 

|l<SU  ■  L  T.  *s- As t. '?</>■%  7 

/I  AT*  -/0-i  X3 


(one/  Subf  i  t le) 

/  'NUMERICAL  GENERATION  OF 
ORTHOGONAL  &JRVILINEAR  c£ 
JUCLIDEAN  ^PACE  .  ^ 


1IMENSI0NAL 
^TES  IN  AN 


7.  authorm; - - _  - - - - _. -  - 

I  Z.  U.  A.  Warsi  R.  A.^Weed  J.  F. /Thompson 


/  0 


9  performing  organization  name  and  aooress 

Mississippi  State  University 
Department  of  Aerospaceengineering 
Mississippi  State,  MS  39762 _ 


II.  CONTROLLING  OFFICE  NAME  AND  ADDRESS 


.UM  I  nWUL'NV  r  IVt  nNmi.  mru  ni/w.wM-  . 

Air  Force  Office  of  Scientific  Research  /A'/ft 
Bolling  AFB,  Washington,  X  20332 


RK AD  INSTRUCTIONS 
UK1-ORK  COMPLETING  hOR.M 
3.  RECIPIENT' 


5.  TYPE  OF  RE 


&  PERIOD  COVERED 


(Q  .  Interim  PtL.pT 

pXqf  fllRMiNG  O^G.  REPORT  NLU 


sat 


8.  contract  or  grant  number^; 


A  AFOSR-^-^185 


7 


'10.  PROGRAM  ELEMENT.  PROJECT.  T  ASK 
AREA  6  WORK  UNIT  NUMBERS 


-  ^  from  Controlling  Ollice) 


kl^nmETTU  r  PAGE  s 

39_ 


IS  SECURITY  CLASS.  (r> I  Ih u  report) 


^UNCLASSIFIED. 


IS*  DECLASSIFICATION  DOWNGRADING 
SCHEDULE 


16.  DISTRIBUTION  STATEMENT  <ot  thts  Report) 


Approved  for  public  release;  distribution  unlimited 


17  DISTRIBUTION  STATEMENT  (ol  the  abstract  entered  it'  fllock  20,  It  lilferent  from  Hopnrf) 


10  SUPPLEMENTARY  notes 


19  KEY  WORDS  (Confmuo  op  reverse  <»>/*’  tf  ne<  e  ***trv  and  itienliiv  hy  block  number) 


Grid  Generation,  Curvilinear  Coordinates,  Numerical  Methods,  Computational 
Fluid  Dynamics 


20  ABSTRACT  (Continue  mt  reverse  aide  It  neceaamry  and  identify  h.  hlo<  fc  numhur) 

•'  In  this  paper  a  non-iterative  method  for  the  numerical  generation  of 
orthogonal  curvilinear  coordinates  for  plane  annular  regions  between  two 
arbitrary  smooth  closed  curves  has  been  developed.  The  basic  generating 
equation  is  the  Gaussian  equation  for  an  Euclidean  space  which  has  been 
solved  analytically.  The  method  has  been  applied  in  many  cases  and  these 
test  results  demonstrate  that  the  proposed  method  can  be  readily  applied 
to  a  wide  variety  of  problems.  ._y 


T 


UNf 


TASSIFIKD  2  ?Q  J.  1 


SECURITY  CL  ASSt  r  1C  ATION  O r  TmS  PAGE  '*t<en  Pat*  kotered) 


OU  X  &  £*  O  V 


DD  ,7.:“,  1473  EDITION  OF  I  NOV  65  15  08SOLETE 


NUMERICAL  GENERATION  OF  TWO-DIMENSIONAL  ORTHOGONAL 
CURVILINEAR  COORDINATES  IN  AN  EUCLIDEAN  SPACE 


.j'J;,  ij  j 


Z.  U.  A.  Warsi,  R.  A.  Weed,  and  J.  F.  Thompson 

Department  of  Aerospace  Engineering 
Mississippi  State  University 
Mississippi  State,  MS  39762 

Summary 

In  this  paper  a  non-iterative  method  for  the  numerical  generation 
of  orthogonal  curvilinear  coordinates  for  plane  annular  regions  between 
two  arbitrary  smooth  closed  curves  has  been  developed.  The  basic  generating 
equation  is  the  Gaussian  equation  for  an  Euclidean  space  which  has  been 
solved  analytically.  The  method  has  been  applied  in  many  cases  and  these 
test  results  demonstrate  that  the  proposed  method  can  be  readily  applied 
to  a  wide  variety  of  problems. 


Accession  F 

NT IS  GRAAI 
DTIC  TAB 
Unannounced 
Justif icati< 

»r 

□ 

□ 

5n_ 

By 

Distribution/ 

I  Avai 

iDist 

ft 

labilil 

Avail 

Spec 

.y  Codes 

and/or 

ial 

$ 


DTIC 

ELECTE 
OEC  30  1980 


D 


AiK  FORCE  OFFICE  OF  SCIENTIFIC  RESEARCH  (AFSC; 
NOTICE  OF  TRANSMITTAL  TO  DDC 
This  technical  report  has  been  reviewed  and  is 
approved  for  public  release  I AW  APR  190-12  (7b) • 
Distribution  is  unlimited. 

A.  D.  BL0S1 


Technical  Info 


"gin  2  29  041 


J  ■  a 


Introduction 


The  problem  of  generating  orthogonal  or  non-orthogonal  curvilinear 
coordinate  systems  in  arbitrary  domains  is  a  problem  of  current  interest 
in  many  branches  of  physics  and  engineering,  and  particularly  in  fluid 
mechanics  and  aerodynamics.  The  idea  of  generating  coordinate  meshes 
by  numerically  solving  a  set  of  partial  differential  equations  under  the 

boundary-geometric  data  as  the  boundary  conditions  arose  with  the  work  of 

1  2  3  A  3 

Winslow  .  Later  Barfield  ,  Chu  ,  Godunov  and  Prokopov  ,  Amsden  and  Hirt 

and  Potter  and  Tuttle**  used  this  concept  in  generating  coordinate  systems 

for  particular  physical  situations.  The  whole  concept  has  however  been 

used  in  a  much  more  organized  manner  by  Thompson,  Thames  and  Mastin^  (later 

Q 

referred  to  as  the  TTM  method)  in  developing  and  coding  the  computer 
program  for  generating  non-orthogonal  coordinates  in  a  variety  of  two- 
dimensional  situations.  The  user,  however,  has  no  control  over  the 
orthogonality  or  non-orthogonality  of  the  generated  coordinates. 

The  underlying  basis  of  all  the  above  methods,  including  that  of 

9  10  11  12 

Pope  ,  Starius  ,  Middlecoff  and  Thomas  ,  and  Mobley  and  Stewart  is  the 

choice  of  a  set  of  coupled  partial  differential  equations.  Two  exceptions 

13 

to  the  above  are  the  methods  of  Eiseman  ,  whose  method  is  of  an  algebraic- 

14 

geometric  nature,  and  that  of  Davis  which  is  based  on  the  Schwarz- 
Christoffel  transformation  of  the  complex  variable  theory. 

In  the  differential  equations  method,  except  for  the  work  of  Starius^ 
where  a  hyperbolic  system  of  equations  is  used,  all  other  methods  rest  on 
the  system  of  elliptic  partial  differential  equations.  The  elliptic 
system  is  usually  a  set  of  Laplace  or  Poisson  equations  V2£  =  -f  and 
V2n  =  -f 2 ,  where  £  =  const,  and  n  =  const,  are  the  coordinate  curves  with 
n  =  const,  on  all  the  boundary  curves.  Looking  a  little  deeper,  one  finds 


that  these  equations  provide  a  set  of  differential  constraints  or  relations 
among  the  fundamental  metric  coefficients  g^,  g^  and  g The  next  step 
is  to  interchange  the  role  of  dependent  and  independent  variables  and  then 
to  solve  the  coupled  system  for  the  Cartesian  coordinates  x  and  y.  In  the 
TTM  method,  the  arbitrariness  of  f^(£,n)  and  f2(£>n)  has  been  used  to  control 
or  redistribute  the  coordinate  lines  in  the  desired  regions. 

In  this  paper  we  develop  a  new  approach  based  on  providing  another 
relationship  among  the  fundamental  metric  coefficients  which  is  not  based 
on  any  arbitrary  assumption.  This  relationship  is  provided  by  the  condition 
that  the  coordinates  are  to  be  generated  in  an  Euclidean  space.  The  most 
natural  choice  is  then  to  use  the  Gaussian  equation^"*  for  an  Euclidean 
space,  viz.,  a  space  of  zero  curvature.  This  fundamental  equation  is  one 
equation  in  the  three  unknowns  g^,  g^  and  To  close  this  equation 

we  can  use  the  simplest  elliptic  system  of  two  Poisson  equations. 

The  preceding  ideas  have  been  tested  in  the  generation  of  orthogonal 
coordinates  in  the  annular  region  between  two  arbitrary  smooth  closed 
curves.  In  the  case  of  orthogonal  coordinates  (i.e.,  g^?  =  0)  the  resulting 
equations  show  that  g£2  is  a  function  of  C,n  and  g^  so  that  a  single 
equation  for  the  determination  of  =  x^  +  y^  is  obtained.  A  general 
solution  of  this  equation  with  g^  =  g^2  can  be  written  down  in  a  series  form 
with  the  Fourier-coef f icients  determined  from  the  prescribed  values  of  g^ 
at  the  inner  and  outer  boundaries.  Further,  from  the  earlier  work  of  Potter 
and  Tuttle^  we  have  the  result  that  in  the  case  of  orthogonal  coordinates 
the  ratio  8^^822  3  Pro^VKt  functions  of  i.  and  n.  This  result  can 


be  used  to  devise  new  coordinates  .f, '  and  a '  in  which  the  resulting  equation 


is  again  of  the  same  form  as  in  f,  and  n.  Thus  the  same  solution  can  be  used 


with  a  change  of  variables  to  provide  the  solution  when  g^2  ^  either 

with  or  without  coordinate  redistribution. 

The  method  developed  on  the  preceding  ideas  therefore  provides  a  non¬ 
iterative  closed  form  analytic  solution  for  the  case  of  two-dimensional 
orthogonal  coordinates.  The  main  methodology  is  detailed  in  the  succeeding 
sections.  Numerical  results  of  the  generated  coordinates  are  shown  in 
Figures  3-8. 


I 


Formulation  of  the  Problem 


For  the  purpose  of  continuity  of  presentation,  we  first  state  the 
following  well  known  results:  the  line  element  ds  in  any  space  is  given 
by  the  Riemannian  formula 


(ds)2  =  gijdxidx-'  (1) 

where  the  g  's  are  c^le  covariant  components  of  the  fundamental  metric 

k 

tensor.  The  choice  of  the  coordinate  system  (x  )  for  any  space  is  quite 
arbitrary,  viz.,  any  coordinate  system  can  be  introduced  for  reference 
purposes,  however,  the  values  of  g (x  )  and  their  distributions  "in  the 
small"  depend  both  on  the  coordinate  system  and  on  the  intrinsic  geometry 
of  the  space  in  such  a  manner  that  the  same  value  of  ds  is  obtained 
irrespective  of  the  chosen  coordinate  system.  Further,  if  in  the  chosen 
space  it  is  possible  to  introduce  a  set  of  rectangular  coordinate  axes, 
then  the  line  element  is  also  given  by 


(ds)2  =  6..dx1dxj  (2) 

ij 

k 

where  is  the  Kronecker  delta  and  (x  )  is  now  a  rectangular  Cartesian 
ij 

system.  A  space  in  which,  in  addition  to  any  general  coordinate  system,  the 
line  element  is  also  obtainable  through  equation  (2),  is  known  as  an 
Euclidean  space. 

The  preceding  ideas  can  be  condensed  by  introducing  the  concept  of 
curvature  of  the  chosen  space  in  which  the  coordinate  system  (x  )  has  been 
introduced.  If  the  Riemannian  curvature  of  a  space  is  zero  then  the  space 
is  said  to  be  Euclidean  and  it  is  then  possible  to  introduce  rectangular 
Cartesian  coordinates  in  this  space.  Therefore,  if  a  two-dimensional  space 


is  Euclidean  then  t lie  Riemannian  or  Gaussian  formula  expressing  the  relation 


among  the  g^.'s  f°r  any  coordinate  system,  (while  writing  x1  =  4,  x2 
is  given  by 


n) , 


X  .  X  .  0 


3n  g 
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34  '  g 


(3) 
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where  r  are  the  Christoffel  symbols  of  the  second  kind  defined  by 


i  ,  i  n 

lk  2  V  9x1  9x‘ 


and 


811  822  ~  (g12)  * 


(4) 


The  main  theme  of  the  present  paper  is  the  choice  of  the  fundamental 
equation  (3)  for  the  determination  of  g„'s  and  then  to  determine  the 
rectangular  Cartesian  coordinates  x  and  y  as  functions  of  4  and  q. 


Determination  of  x  and  y 

Equation  (3)  implies  that  there  exists  a  continuous  function  a(4,n) 
such  that 


a{  =  ^  r2 

4  gu  11 


a  =  r2 
n  gu  12 


(5) 


where  a  variable  subscript  denotes  partial  differentiation.  Based  on  the 
following  formulae 


811  =  xl  +  yl'  812  =  Vn  +  Vn’  822  =  Xn  +  yn 


ril  =  7=  (V44  "  y4X44}>  f12  =  7=  (x4y4n  '  Vm 

*e  *  g 


(6) 
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I 
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it  is  easy  to  show  by  direct  substitution  that 


cos  a,  y ^  =  -✓g11  sin  a 


x  =  — —  (g.  „  cos  a  +  /g  sin  a)  ,  y  =  — - —  (/g  cos  a  -  gn  sin  a)  /  (7) 

n  ZZ  12  n  4  12 

gll  gll 


Using  equations  (5)  and  (7)  in  the  expressions  for  da,  dx  and  dy,  we  obtain 
the  expressions  for  a,  x  and  y  (first  obtained  by  Martin1^). 


*  "  ~  /  (r2ll  +  r212  d"> 


(8) 


x  =  / l  '''g  1 1  cos  a  d(  +  ~rr —  (g  „  cos  &  +  ^  sin  a)dn] 

11  /&H  LZ 


(9) 


y  =  / f_v/8n  sin  a  d^  +  ~7^~  cos  a  -  St  9  sin  °0drl] 


Zg- 


11 


12 


(10) 


The  geometrical  interpretation  of  a  is  that  it  is  the  angle  of  in¬ 
clination  with  respect  to  the  x-axis  of  the  tangent  to  the  coordinate  line 
n  =  const,  directed  in  the  sense  of  increasing  values  of  the  parameter  £,• 

The  choice  of  the  minus  sign  in  equation  (5)  is  due  to  the  adopted  convention 
that  both  4  and  a  be  treated  as  positive  in  the  clockwise  sense. 


Case  of  Orthogonal  Coordinates 


Equation  (3)  when  written  in  terms  of  g_  and  g  =  g^g22  ~  (g^)' 


has  the  form 


812  3gll 


1  3g22 


34 


gj^'g  3n 


,  J_  ,2  3g12 

3n  '  r  . 


1  3gll 


"g  3f, 

g12  3gll 


)  =  0 


(ID 


g  34  ✓  g  ,)!) 


gxl  ►g  34 
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Here  a 


which  is  an  equation  in  three  unknowns,  viz.,  g  ,  g  and 
wide  range  of  possibilities  are  open  to  express  g^  and  g ^  as  functions 
of  g^  either  through  algebraic  or  differential  relations.  Fortunately, 
in  the  case  of  orthogonal  coordinates  this  arbitrariness  is  minimal  and 
equation  (11)  can  be  reduced  to  a  very  simple  form  which  can  be  solved 
analytically.  The  following  discussion  and  analysis  pertains  only  to 
orthogonal  coordinates. 

In  the  case  of  orthogonal  coordinates,  the  coefficient  g^  is  zero, 

i.e.  , 

g12  =  Vn  +  y€yn  =  0  (12) 

which  is  satisfied  by  the  equations 


where  F(x,y , £ , n ,x. ,y^ ,x^ ,y^)  >  0  is  a  continuous  function  of  its  arguments. 

Let  the  boundary  of  a  bounded  region  ft  in  an  Euclidean  two- 
dimensional  space  be  a  simple  smooth  curve  x  =  xm(£)  ,  y  =  y^  (£)  >  with  a 
uniformly  turning  tangent.  In  the  region  ft,  let  ftg  be  an  annular  subregion 
bounded  by  the  inner  boundary  1’^  and  the  outer  boundary  T  as  shown  in 
Figure  1.  The  region  ftg  is  to  be  mapped  onto  a  rectangular  region  R  in  the 
tin-plane  by  a  transformation  with  the  initial  data 


I 


so  as  to  have 


x  =  *JO 

y(f.»n J  =  yM  ( 4 ) 

x  =  x(;„,n) 

y  =  y(i,n)  ' 


0  <  £  <  £ 

—  m 


(14) 


n£l  -  n  --  "o 


(15) 


•  v  •('*  * .  V**  » 


:■  ?3*r  w 


where  n  and  n  are  the  actual  parametric  values  associated  with  the 

B  “ 

boundaries  S'^  and  ^ ,  respectively,  and  x  and  y  are  periodic  in  the  £- 
argument  with  the  period 

21  =  £  (16) 
m 

The  set  of  equations  (13)  along  with  the  initial  data  (14)  can  form  a 
well-posed  initial-value  problem  for  hyperbolic  equations  if  certain 
conditions  on  F  are  satisfied.  It  is  to  be  noted  that  in  this  case  no 
boundary  data  is  needed  on  the  curve  1’^.  This  problem  has  been  considered 
by  Starius-*-^ 

Since  in  this  paper  the  equation  for  the  condition  of  an  Euclidean 
space  (equation  (11))  forms  the  basis  of  the  proposed  method,  we  expect 
to  have  an  elliptic  boundary-value  problem  to  be  solved  under  the 
Dirichlet  conditions.  Equation  (11)  with  the  substitution  g^  =  0  and 
on  using  equations  (13)  becomes 


_3_ 

3 1, 


Fg 


i - —  (F2g  )]  +  .'L  [_I 

3^  (  girJ  3n  Fg 


3g 


11, 


11 


11 


3n 


0 


(17) 


8?2  ^  gll 

g  =  (f811)'  (IS) 

Equation  (17)  can  now  be  solved  as  an  elliptic  boundary-v3lue  problem 
in  g^  provided  that  it  can  be  proved  that  F  is  a  function  of  4,n  and  8-q- 
These  considerations  on  F  will  also  provide  a  class  of  functions  from  which 
F  can  be  chosen  in  a  simple  way. 

In  Reference  10,  Starius  has  proved  the  following  properties  of  F: 

(a)  F  is  an  invariant  under  translation  and  rotation  of  the  coordinate 
axes  (x,y).  Therefore,  if  (x,y)  is  a  solution  of  equations  (13),  then  it 


9 


can  be  shown  that 


x  =  a  +  x  cos  6  -  y  sin  6 
y  =  b  +  x  sin  6  +  y  cos  6 
is  a  solution  of  the  equations 


y  =  F  x 
n 


where  a,b  and  0  are  arbitrary  constants.  These  considerations  show  that  F 
is  not  an  explicit  function  of  x  and  y. 

(b)  On  the  basis  of  the  results  obtained  in  (a) ,  we  have 


F(C,n, 


F(C,n,  xc,  yv  xn,  yn) 


hence  —  =  0  for  all  values  of  0  including  0=0.  Evaluating 
clU 

(4fo  ,  we  obtain 
d0  0  =  0 


9F  ,  3F  9F  3F 

3x_  XC  3y„  3x  Xq  3y 

4  i  n  n 


o 


This  equation  shows  that  F  =  F(£,n,  S^l’  §22^’  so  t^at  F  depends  on 
Xjr  >  y  y^  through  gn  and  z21  only. 


(c)  The  considerations  in  (a)  and  (b)  along  with  equation  (18)  show  that 


§22  ^11  ^11*  §22^ 

so  that  in  principle  g^  can  be  expressed  as  a  function  of  g^  and  con¬ 
sequently  F  =  F(f,,n,  g^) .  This  proves  the  contention  that  in  the  case 
of  orthogonal  coordinates,  equation  (17)  is  sufficient  for  the  calculation 

°f  gir 

The  set  of  functions  F  >  0  satisfying  the  conditions  (a),  (b) ,  and 
(c)  enumerated  above  also  contain  F  =  1  as  an  element.  This  choice  of  F' 


*This  is  by  no  means  a  restriction,  as  is  demonstrated  later. 


yields  the  simplest  possible  form  of  the  generating  equation.  With  F  =  1, 
equation  (17)  becomes 


where 


£p_  +  .  o 

H2  3n2 


P  =  in  gu,  g22  =  g1;L,  V2i  =  0** 


The  boundary  conditions  are 


(19) 


(20) 


P  =  P  (O  at  n  =  0 


=  P  (£)  at  n  =  n 


>0  <  C  < 


(21) 


where  the  subscripts  3  and  "  denote  the  inner  and  outer  boundaries 
respectively.  The  periodicity  requirement  is  that 


P(C,n)  =  PU  +  24,  n)  (22) 


where 


2i  =  4  • 

m 

A  general  analytic  solution  of  equation  (19)  under  the  conditions  (21) 
and  (22)  is 

-  r  Win 

—  00  nTi  n  tt  f  nTT  l  n  .  00 

P(4,n)  =  a  +  nK  +  X-  sinh  -5-  (n  -n)(a  cos  — —  +  b  sin  — —  )/sinh  — 5— 
o  n  j.  n  Jc  n  x,  k. 


n7in 

+  t  sinh  —21  (c  cos  j  sin  )/sinh  — — 

n=i  S.  n  S.  n  8.  <■ 


(23) 


whe  re 


K  =  (c  -a  )/n 
o  o  “ 


(24) 


PThere  is  no  loss  of  generality  in  setting  the  parametric  value  u  .  =  0.  The 
value  n  must  be  interpreted  as  the  difference  between  the  actual  values  at 
the  outer  and  inner  boundaries.  The  determination  of  nm  is  of  crucial 
importance  to  this  work  and  is  discussed  in  the  next  section. 

!<*Refer  to  the  next  section. 


1  1  21 
*o  =  U  f o  ca  =  u  fQ  Pa,(t)d£. 


*n  =  I  C  V° 


d^*  bn  =  K  V«>  ^^dc}(25) 


,  21  e  21  1 
cn  =  I  /Q  P„(0  cos  dCf  dn  =  -  !o  ?JO  sin  dc  I 

Having  determined  the  coefficients  a  ,  b  ,  c  and  d  as  defined  in  (25), 
°  n  n  n  n 

we  can  obtain  the  values  of  g^  from  (23)  for  all  values  of  £  and  q.  To 
find  the  expression  for  a  we  consider  equations  (5),  which  for  orthogonal 
coordinates  are 


1  9sll  1  9g22 

=  2/g  “n  =  "  2/i  « 

For  g22  =  g^,  these  equations  become 

1  3P  _  1  3P 

ac  =  2  3n’  ”  2 

and  on  integration  yield  the  exact  expression  for  a. 


a(S,n)  =  a(£,o)  +  T.  — — — - 

n=l  nun 

2  sinh  — 


cosh  (n  -n)  , 

*.  00  /t  nn£  ■_  mif,, 

-  (bn  cos  T"  "  3n  Sin 


+  T. 

n=l 


2  sinh 

nun 

oo 

l 

MU} 

CD 

cosh  ^ 

2  sinh 

nrrn 

oo 

l 

1 

(c  sin 


nH  nnf,  , 

— —  -  d  cos  — —  ) 


,,  nnf,  .  nnf,. 

(b  cos  — - a  sin  —r—) 

n  2  n  2 


-  (c  sin  — - - d  cos  z  ) 

nun  n  l  n  l 


2  sinh  — • 


Since  the  line  integrals  for  the  determination  of  x  and  y  (cf.  Eqs.  (9) 
and  (10))  are  independent  of  the  path,  viz., 
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!  i 


3  3 

—  (/iu  cos  a)  =  —  (/i22  sin  a) 


r 


hence 


—  (/in  sin  a)  =  -  --  (/i22  cos  a) 


c(C,n)  =  x(C,o)  +  Jq  /g22  sin  a  dr 


(28) 


(29) 


y(S,n)  =  y(£,o)  +  jQ  /g22  cos  a  dn 

The  preceding  analysis  completes  the  basic  development  of  the  subject. 


Coordinate  Re-Distribution  (Contraction) 

In  order  to  have  the  capability  of  re-distributing  the  coordinate  lines 
so  as  to  have  a  control  on  the  mesh  spacings  in  the  desired  regions,  we 
consider  a  transformation  from  (£,n)  to  new  coordinates  (£,n).  The  trans¬ 
formation  functions  can  arbitrarily  be  selected,  but  can  also  be  linked  in 
some  manner  to  the  physical  field  behaviors.  For  example,  in  viscous 
flow  problems  the  effect  of  viscosity  near  a  wall  can  be  incorporated  in 
the  transformation  functions.  Below  we  proceed  without  specifying  these 
functions  and  then  give  one  example  in  the  section  on  numerical  method. 

On  transformation  from  (£,,n)  to  (i,n),  the  covariant  metric  coefficients 
transform  as 


8U 


g 


,  k  a  l 

dX  dX 

kH  -,-i  T~j 
3x  3x 


so  that  on  using  the  relation  g22  =  g^  and  g^2  =  0,  we  have 


!i 


We  now  introduce  the  transformation 


4  =  HU  I 

n  *  f (n)  >  (31 ] 
where  the  functions  <p  and  f  are  continuously  differentiable  and  satisfy  the 
condit ions 


*(4  )  -  o,  f(n  )  «  0 

O  fcs 


where  4  *  0  and  n  *  0  correspond  respectively  to  4 

v  =  di  ,  „  =  if. 
dt.  dn 


=  4  and  n  ■  n  •  Defining 

o  t- 


Ula) 


we  have  from  (30) 


g1l (4 . n)  =  ■'  g11(4,n) 


g22(4,n)  =  ^t‘g11(4,ri) 


i, , (4,n) 


(32b) 


( 32c ) 


(32d) 


To  obtain  the  solution  in  the  (1,7))  coordinate  system,  we  merely  have 
to  replace  f,  and  n  by  the  functions  ip(4)  and  <Kn)  respectively  in  (23)  and 
(27),  while  (28)  and  (29)  become 


X (4 , n)  =  X (4 , f.Q)  +  /  '  ••  g22  sin  a(4,n)dn 


Y(4,n)  =  Y(4,4Q)  +  Jn  <t22  cos  ll(4 , n)dii 


(33a) 


(33b) 


it  must  be  noted  that  on  transformation  the  resulting  metric  coefficients 
and  g22  are  not  equal. 

The  salient  feature  of  the  preceding  analysis  is  that  the  solution 
under  the  condition  g^  =  g^  can  be  used  to  obtain  the  solution  for  the 
case  g^2  +  by  coordinate  transformation. 


» ***  -A 


Uniqueness  Condition  on  t,  for  Orthogonal  Coordinates 

Any  method  for  the  generation  of  orthogonal  coordinates  on  the  preceding 
lines  has  to  be  supplemented  with  a  uniqueness  condition  on  the  behavior 
of  4  and  a  method  for  its  selection.  The  following  analysis,  besides 
covering  the  above  two  aspects,  also  provides  a  general  basis  for  the 
earlier  choice  F  =  1. 

For  the  case  of  orthogonal  coordinates,  equations  (7)  can  also  be 
expressed  In  the  inverse  form  as 


•  x  »  cos  ,  f.y  *  -sin 


hx  =  sin  a/.g22  ,  »y  =  cos  u/.'g 


Writing  F  =  •  g.r,7"gj'j ,  eliminating  i  between  the  above  equations  and  then 
bv  cross  differentiation,  we  obtain  the  following  equations: 


(  */F)  +  ,v  VF)  =  ° 


>x  X 

(F"X> + <Fv =  ° 


(34) 

(35) 


Following  Trotter  and  Tuttle*’  we  assume  that  the  (.-curves  in  the  xv-planc 
are  free  1 rom  sources  and  sinks.  This  condition  establishes  a  unique 
correspondence  between  the*  points  on  each  pair  of  u  —  constant 
lines.  In  the  absence  of  sources  or  sinks,  we  have 


di  v|  grad  ,*(■»)!  =  0 

where  ,(*i)  is  an  arbitrary  differentiable  function  of  n,  and  as  such  grad 
;,(,,)  is  oriented  along  the  normal  to  the  curve  n  =  const.  Using  the 
expressions 

i  grad  (i,  »  1  / ►  fi 2 .  V- n 
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in  (3b) ,  we  obtain 


—  (in  /r~/T“)  m.  /  *i- 

an  ^  n  gll/g22;  dn2  dn 


Writing  ~  =  l/v(n)  and  denoting  the  arbitrary  function  due  to  integration 


as  InuiO »  we  obtain  the  result 


*811^822  “  f  ~  v(n) 

Introducing  new  variables 


(37) 


'  -  /-«>«.  »-  •  /  ^ 


(38) 


and  using  (37)  in  (34)  and  (35),  we  get 

V*  (,  '  =  0  (39) 

V-V  =  0  (40) 

Further  using  (37)  and  (38)  in  the  fundamental  equation  (17),  we  obtain 


=  o 

at.  *•’  an"’ 


(41) 


where 


P‘  =  fn  *11 


8  >'>  ~  S 


11 


gu  =  x-  .  +  y;.  .,  g2,  =  x;r  +  y;. 


The  solution  of  equation  (41)  is  of  the  same  form  as  that  of  equation 
(19),  viz.,  (23),  and  is  obtained  by  replacing  \  and  n  bv  '  and  n‘. 
However,  the  important  result  obtained  here  is  that  a  generating  equation 
of  the  form  (19)  or  (41)  must  be  supplemented  with  a  Laplace  equation  for 
or  respectively.  This  is  the  result  or  the  required  uniqueness 
condition  on  >,  for  the  generation  of  orthogonal  coodinates.  The  condition 
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on  r)',  i.e.,  equation  (AO),  is  implicitly  satisfied  by  the  coupled  system 
of  equations  (39)  and  (41)  and  needs  no  discussion. 

The  condition  equation  V2£  =  0  can  rigorously  be  satisfied  in  all 
cases  if  we  take  e,  as  the  angle  traced  out  in  a  clockwise  sense  by  the 
common  radius  of  the  concentric  circles  in  a  conformal  representation  of 
the  inner  and  outer  boundaries.  The  numerical  scheme  for  this  aspect  of 
the  problem  is  an  iterative  one.  In  place  of  this  elaborate  scheme  we  have 
devised  another  method  which  is  much  simpler  and  non-iterative.  Both 
of  these  methods  have  been  discussed  in  the  next  section. 


Numerical  Method  of  Solution 


Based  on  the  formulation  of  the  problem  as  discussed  in  the  preceeding 
sections,  we  now  have  a  non-iterative  algebraic  computational  problem 
which  can  be  handled  in  a  straight  forward  manner.  However,  before  solving 
a  specific  problem,  it  is  important  first  to  establish  an  orthogonal 
correspondence  between  unique  points  of  the  inner  and  outer  boundary 
curves  which  are  to  be  connected  by  a  specified  number  of  i,  =  constant 
curves,  and  second  to  obtain  the  numerical  value  of  the  parametric 
difference  n,,.  Two  methods  for  the  establishment  of  x((.)  and  y(0 
are  given  below. 

Method  1:  It  has  been  mentioned  in  equation  (20)  and  later  discussed 
in  the  preceeding  section  that  the  curves  =  constant  in  the  physical  xy- 

plane  must  satisfy  the  Laplace  equation  V  =  0.  For  this  condition  to 

be  satisfied  we  can  take  as  the  angle  traced  out  in  the  clockwise  sense 
by  the  common  radius  of  the  concentric  circles  in  a  conformal  representation 
of  the  inner  and  outer  boundary  curves  as  follows. 

The  function  z  =  t  tz  ’),  which  conformally  maps  the  region  of  the  z- 
plane  exterior  to  the  specified  curve  ('  onto  the  region  of  the  z  ’-plane 
exterior  to  t  lie  circle  C'  of  radius  a,  can  be  represented  hv  a  Laurent's 
expansion  as 

z  =  z  '  +  p  +  iq  f  :  (p  +  iq  )l‘\)n  (42) 

o  o  .  u  n  z 

1 

For  points  on  the  circumference  of  the  circle  C 

z  ’  =  a  e 


1H 


.  »  ^  ^ 


so  that 


xU)  "  P0  +  (Pj+a)  cos  4-  sin  4  +  l  (pn  cos  n4  -  q  sin  nC)  (43) 

n=2  n 

y(4)  =  q  +  (p,-a)  sin  4+  q.  cos  4  +  I  (p  sin  n4  +  q  cos  n4)  (44) 

±  „n  n 

n=2 

The  same  form  of  equations  can  be  written  for  the  outer  boundary  with 
A  as  the  radius  of  the  circle  in  the  conformal  plane. 

Now  y  =  y(x)  is  known  either  in  functional  or  tabular  form.  Thus 
starting  from  an  initial  guess  for  x  the  corresponding  ordinates  are  used 
to  determine  the  Fourier  coefficients  from  equation  (44)  which  in  turn 
determine  a  new  set  of  abscissae  and  then  the  ordinates,  and  so  on.  The 
convergence  of  this  iterative  method  yields  x  =  x(4)  and  y  =  y(4)  both 
for  the  inner  and  outer  boundaries.  Note  that  after  the  completion  of 
convergence  we  have 

a  =  2^  lo"  1x^(4)  cos  4  -  y^C)  sin  (45a) 

and 

A  =  /2"  [x  (4)  cos  4  -  y  (4)  sin  4 ]d4  (45b) 

Method  II:  In  lieu  of  using  Method  I,  we  have  obtained  equally  good 

results  by  proceeding  as  follows.  This  method  looks  to  be  equivalent  to  Method  I. 

The  inner  and  outer  boundary  data  is  available  to  us  either  in  tabular 
or  functional  form  as 

yp  =  y(xtJ,  y„  =  y(xj  (46) 

We  now  circumscribe  circles  around  the  inner  and  outer  boundary  curves.  Two 

cases  arise  depending  on  whether  the  circles  are  concentric  or  nonconcentric . 

Case  I:  If  the  circumscribed  circles  are  concentric  (Fig.  2a),  then 

we  select  those  sets  of  ordinates  which  correspond  to  the  abscissae 
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Determination  of  n 

oo 

The  parametric  difference  r)  is  connected  in  some  manner  with  the 

CO 

"modulus"  of  the  domain.  Determination  of  the  modulus  for  annular  regions 

has  been  considered  by  Burbea^  and  Gaier^.  In  this  paper  we  base  the 

determination  of  on  the  radii  a  and  A  in  the  conformal  representation 

process  described  in  equations  (45),  and  define  it  as 

n  =  fn(^)  (48) 

®  'a 

If  Method  I  is  used  then  A  and  a  are  available  as  parts  of  the  converged 

iteration,  while  if  Method  II  is  used,  then  they  are  obtained  by  simple 

quadratures  applied  to  equation  (45). 

Having  determined  the  appropriate  x(£,)  and  y(£)  both  for  the  inner 

and  outer  boundaries,  we  first  calculate  the  values  of  (g.  )  and  (g  ) 

i.  1  P  11  w 

numerically  and  then  of  P„(£)  and  P  (O .  Based  on  these  distributions 

P  00 

the  Fourier  coefficients  a  ,  b  ,  c  and  d  are  computed  by  numerical 
quadrature  through  the  use  of  equation  (25).  Since  these  values  of  the 
coefficients  are  independent  of  the  spacings  between  n  =  constant  lines 
the  same  values  are  used  when  a  redistribution  of  n  =  constant  lines  is 
desired.  Substituting  the  Fourier  coefficients  in  equation  (23)  we 
determine  the  values  of  P((,,n)  and  hence  of  g^(£,n)  for  the  whole 
annular  region.  A  knowledge  of  P(t,,n)  determines  a(f,,n)  through  equation 
(27)  and  so  also  the  Cartesian  coordinates  x(£,,n)  and  y(£,n)  through 
equations  (28)  and  (29)  by  numerical  quadrature. 

A  computer  program  with  the  option  of  redistributing  the  coordinate 
lines  in  any  desired  manner  lias  been  written  and  used  to  generate  the 
orthogonal  coordinates  for  various  annular  regions  as  shown  in  Figures  3-8. 
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For  the  example  problems  we  have  chosen  the  following  forms  of  the  function 


ip  and  f  of  equations  (31) . 


'(C)  = 


C  -c 

m  o 


so  that 


C  ~C 

m  o 


hc(n-nft)  K^n_rV 
f  (n)  =  ^ - 

n“"nB  A'V 


-  -  K(n_nR) 

- — —  [l+(n-n  )«.nK]  — - e_ 

n»~n8  K(^~V 


where  K  >  1  is  an  arbitrary  constant,  £  =  n ,  and  C  =  C^,  n  =  correspond 

respectively  to  C  =  21  =  2tt  and  n  =  n  .  We  treat  C  and  n  as  integers  so 

that  C  =  1,  C  =  IMAX,  n  =  1,  n  =  JMAX.  Since  n  is  known  from  (48), 
o  m  8  00  " 

hence  by  specifying  the  numerical  values  to  K  and  JMAX  we  can  create  the 

desired  mesh  control  in  the  direction  of  n.  The  value  of  K  between  1.05 

20 

and  1.1  is  quite  sufficient  to  have  very  fine  grid  spacing  near  the 
inner  boundary. 

The  number  of  terms  to  be  retained  in  the  series  (23)  is  usually  small 

for  convex  inner  and  outer  boundary  curves,  though  we  have  retained  (IMAX-l)/2 

number  of  coefficients  in  each  computation.  This  number  is  the  optimum 

21 

number  of  terms  in  a  discrete  Fourier  series  having  IMAX  number  of  points 
in  one  period.  The  average  computer  time  for  the  complete  computation  on 
the  UNIVAC  1100/80  for  IMAX  =  73  and  JMAX  =  60  field  is  about  2.75  minutes. 


Summary  of  Numerical  Experimentations 

In  the  course  of  this  investigation  a  number  of  cases  of  inner  and 
outer  boundary  shapes  and  orientations  have  been  tested  through  the 
developed  computer  program.  The  main  conclusions  are  listed  below: 

(i)  The  method  works  very  effectively  for  smooth  and  convex 

boundaries  of  any  shape  and  orientation  (cf.  Figures  ( 3 )  —  ( 5 )  and 
(7)). 

13 

(ii)  For  concave  boundaries  a  method  similar  to  that  of  Eiseman 

has  to  be  used  in  the  placement  of  the  outer  boundary  to  avoid 
intersecting  £  =  constant  lines  (cf.  Fig.  8). 

(iii)  Sharp  turns  or  corners  are  not  admissible  and  have  to  be  rounded 
to  avoid  singularities  (cf.  Fig.  6). 
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Comparison  with  the  TTM  Method 


As  discussed  in  the  introduction,  the  choice  of  the  proposed  method 
has  been  motivated  by  two  considerations:  (i)  to  choose  a  set  of  equations 
which  are  fundamental  to  the  intrinsic  nature  of  the  space  in  which  the 
coordinates  are  sought,  and  (ii)  to  minimize  the  inherent  arbitrariness 
in  the  selection  of  relations  among  the  metric  coefficients.  It  is  the 
purpose  of  this  section  to  show  that  for  the  case  of  orthogonal  coordinates 
the  equations  of  TTM  (Ref.  7,20)  can  be  reduced  in  a  simple  form  which  can 
also  be  used  to  achieve  the  above  two  stated  purposes. 

The  two-dimensional  Laplacian  of  a  scalar  ip  of  general  coordinates  £  and 


n  is 


±  [£  (£22^£ &12 T n 

^  3C  4 


81 -A 


)  (-1-1--n~--12  S] 

3n  ro 


(51) 


Introducing  the  operator 


o2  82235£  2gl23£n  +  8ll3nn 


(52) 


in  (5i)  gives  another  form 

V2 il>  =  -  D2ii  +  ip  V2£  +  ib  V2 
g  a  n  n 


(53) 


i 


Writing  ip  =  x  and  then  =  y  in  (53),  we  recover  the  equation  of  TTM, 
which  are 


D2x  =  -g(xrV2£  +  x^V2n) 


d  y  =  -g(yj-v^£  +  ynv  n) 


For  orthogonal  coordinates  =  0,  so  that 


g22x££  +  gllxnn  =  "r(x£V^  +  \V'n) 


g22y££  +  gllynn  =  'gV?':  +  ynV2n) 
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(54) 


(55) 


ill 


Also,  for  orthogonal  coordinates,  equation  (51)  gives 
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we  obtain 


x.  ,  +  x  .  .  =  0 

4  4  n  n 

y. ..*+y  -  .  =  0 

4  4  n  n 


(59) 


where 


V2£' 


0,  V2!)'  =  0 


From  the  foregoing  analysis  we  conclude  that  the  set  of  equations  (59)  are 
quite  general  and  capable  of  generating  orthogonal  coordinates  (as  has 
earlier  been  proposed  by  Pope^),  and  there  is  no  need  to  solve  a  more 
difficult  and  time  consuming  set  of  non-linear  equations  (55)-  Further 
the  method  of  solving  equations  (59)  follows  the  same  patterns  as  that  of 
the  proposed  method  of  this  paper  as  discussed  before. 
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Conclusions 


A  new  method  for  the  generation  of  orthogonal  coordinates  in  two- 
dimensional  regions  has  been  proposed.  The  method  shows  that  much  of  the 
arbitrariness  in  the  choice  of  relations  among  the  metrical  coefficients 
can  be  minimized  by  the  use  of  the  condition  of  an  Euclidean  space.  The 
method  can  also  be  used  for  simply  connected  regions  only  by  obtaining  the 
solution  of  the  linear  equation  (19)  under  the  changed  boundary  conditions. 
Besides,  the  proposed  method  can  also  be  extended  to  three-dimensional 
regions. 

Acknowledgement 

This  research  has  been  supported  in  part  by  the  Air  Force  Office  of 
scientific  Research,  under  Grant 

/p6rc''^>/7  76  -.7777 


References 


1.  A.  J.  Winslow,  "Numerical  Solution  of  the  Quasi-Linear  Equation  in  a  1 

Non-Uniform  Triangular  Mesh",  Journal  of  Computational  Physics,  2, 

149  (1966).  | 

2.  W.  D.  Barfield,  "An  Optimal  Mesh  Generator  for  Lagrangian  Hydrodynamic  ,  j 

Calculations  in  Two  Space  Dimensions",  Journal  of  Computational  Physics, 

6,  417  (1970). 

\{ 

3.  W.  H.  Chu,  "Development  of  a  General  Finite  Difference  Approximation  •  j 

for  a  General  Domain,  Part  I:  Machine  Transformation",  Journal  of 

Computational  Physics,  8,  392  (1971).  1 

4.  S.  K.  Godunov  and  G.  P.  Prokopov,  "The  Use  of  Moving  Meshes  in  Gas 

Dynamics  Computations",  USSR  Computational  Mathematics  and  Mathematical  ,  J 

Physics,  12,  182  (1972).  I, 

\ 

5.  A.  A.  Amsden  and  C.  W.  Hirt,  "A  Simple  Scheme  for  Generating  General  , 

Curvilinear  Grids",  Journal  of  Computational  Physics,  11,  348  (1973). 

f  i 

6.  D.  E.  Potter  and  G.  H.  Tuttle,  "The  Construction  of  Discrete  Orthogonal 
Coordinates",  Journal  of  Computational  Physics,  13,  483  (1973). 

7.  J.  F.  Thompson,  F.  C.  Thames,  and  C.  W,  Mastin,  "Automatic  Numerical 
Generation  of  Body-Fitted  Curvilinear  Coordinate  System  for  Field 
Containing  any  Number  of  Arbitrary  Two-Dimensional  Bodies",  Journal  of 
Computational  Physics,  15,  299  (1974). 

8.  J.  F.  Thompson,  F.  C.  Thames,  and  C.  W.  Mastin,  "'TOMCAT’-A  Code  for 
Numerical  Generation  of  Boundary-Fitted  Curvilinear  Coordinate  System 
on  Fields  Containing  any  Number  of  Arbitrary  Two-Dimensional  Bodies", 

Journal  of  Computational  Physics,  24,  274  (1977). 

9.  S.  B.  Pope,  "The  Calculation  of  Turbulent  Recirculating  Flows  in 
General  Orthogonal  Coordinates",  Journal  of  Computational  Physics, 

26,  197  (1978).  : 

10.  G.  Starius,  "Constructing  Orthogonal  Curvilinear  Meshes  by  Solving  ; 

Initial-Value  Problem",  Numerische  Mathematik,  28,  25  (1977). 

11.  J.  F.  Middlecoff  and  P.  D.  Thomas,  "Direct  Control  of  the  Grid  Point 
Distribution  in  Meshes  Generated  by  Elliptic  Equations",  AIAA 
Computational  Fluid  Dynamics  Conference,  Paper  No.  79-1462  (1979). 

12.  C.  D.  Mobley  and  R.  J.  Stewart,  "On  the  Numerical  Generation  of 
Boundary-Fitted  Orthogonal  Curvilinear  Coordinate  Systems",  Journal 
of  Computational  Physics,  34,  124  (1980). 


28  fl 


13.  P.  R.  Eiseman,  "A  Coordinate  System  for  a  Viscous  Transonic  Cascade 
Analysis",  Journal  of  Computational  Physics,  26,  307  (1978). 

14.  R.  T.  Davis,  "Numerical  Methods  for  Coordinate  Generation  Based  on 
Schwarz-Christoffel  Transformation",  AIAA  Computational  Fluid  Dynamics 
Conference,  Paper  No.  79-1463  (1979). 

15.  L.  P.  Eisenhart,  "Riemannian  Geometry",  Princeton  University  Press 
(1926). 

16.  M.  H.  Martin,  "The  Flow  of  a  Viscous  Fluid,  I",  Archiv  of  Rational 
Mechanics  and  Analysis,  41,  266  (1971). 

17.  H.  Kofcer,  "Dictionary  of  Conformal  Representations”,  Dover  Publications, 
Inc.  (1952). 

18.  J.  Burbea,  "A  Numerical  Determination  of  the  Modulus  of  Doubly 
Connected  Domains  by  Using  the  Bergman  Curvature",  Mathematics 
of  Computation,  25,  743  (1971). 

19.  D.  Gaier,  "Determination  of  Conformal  Modulus  of  Ring  Domains  and 
Quadrilaterals",  In-Lecture  Notes  in  Mathematics,  No.  399,  Springer- 
Verlag,  Berlin  (1974). 

20.  Z.  U.  A.  Wars!  and  J.  F.  Thompson,  "Machine  Solutions  of  Partial 
Differential  Equations  in  the  Numerically  Generated  Coordinate 
Systems",  Engineering  and  Experimental  Research  Station,  Mississippi 
State  University,  Rep.  No.  MSSU-EIRS-ASE-77-1,  August  (1976). 

21.  F.  B.  Hildebrand,  "Introduction  to  Numerical  Analysis",  McGraw-Hill 
Book  Company,  Inc.,  (1956). 


29 


\\^ 


y  c2 


Cx 


Figure  2.  (a)  Concentric  circumscribed  circles  and  C2  of  radii  rs,  r 

respectively  with  center  at  the  origin.  (b)  Non-concentric 
circumscribed  circles  Cj  and  of  radii  rg  and  r^  and  centers 
at  z  and  z,  respectively.  ^ 


Figure  3.  Confocal  ellipses.  Semi-major  axes  1.48,  5.0,  and  semi-minor 
axes  0.5,  4.802  respectively.  Only  38  n  -  const,  lines  shown 
for  detail. 
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Figure  3.  Confocal  ellipses.  Semi-major  axes  1.48,  5.0,  and  semi-minor 
axes  0.5,  4.802  respectively.  Only  38  n  *=  const,  lines  shown 


Figure  4.  A  blunt  body  section  with  elliptical  outer  boundary. 


Figure  6-  Joukowsky's  airfoil  with  slightly  rounded  trailing  edge. 


Figure  7.  Non-concentr ic  ellipses.  Size  data  same  as  in  Figure  3. 

z  =  (0,0),  z  (1,0). 
s  L 
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